Critical entropies for magnetic ordering in bosonic mixtures on a lattice 
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We perform a numeric study (worm algorithm Monte Carlo simulations) of ultracold two- 
component bosons in two- and three-dimensional optical lattices. At strong enough interactions 
and low enough temperatures the system features magnetic ordering. We compute critical temper- 
atures and entropies for the disappearance of the Ising antiferromagnetic and the icy-ferromagnetic 
order and find that the largest possible entropies per particle are ~ 0.5fcs. We also estimate (op- 
timistically) the experimental hold times required to reach equilibrium magnetic states to be on a 
scale of seconds. Low critical entropies and long hold times render the experimental observations of 
magnetic phases challenging and call for increased control over heating sources. 

PACS numbers: 67.85.Hj, 67.85.Fg, 67.85.-d 



I. INTRODUCTION 

At the moment, one of the prominent focuses and ma- 
jor challanges of experiments with ultracold gases is the 
realization of configurations which can be used to study 
quantum magnetism [l|, Q ■ Though interesting and fun- 
damental on its own, better understanding of (frustrated) 
magnetic systems is further motivated by its relevance 
to high-T c superconductivity and applications to quan- 
tum information processing. Direct studies of condensed 
matter spin systems experimentally are limited by the 
lack of control over interactions, geometry, frustration, 
and contaminating effects of other degrees of freedom. A 
new approach consists of using ultracold atoms in op- 
tical lattices (OL) provided that the system is driven 
towards regimes where it is possible to map the corre- 
sponing (Bose-)Hubbard Hamiltonian to spin models. 

Striking advances in experimental techniques, e.g. high 
controllability and tunability of Hamiltonian parameters, 
and, more recently, single site and single particle imaging 
@j m IE HI ! brought forward the idea, originally proposed 
by Feynmann, of quantum simulation/emulation [7[. In 
the last decade, a considerable amount of theoretical and 
experimental research has been devoted to the objective 
of using ultracold lattice bosons and fermions to address 
many outstanding condensed matter problems via Hamil- 
tonian modeling. Perhaps the biggest remaining exper- 
imental challenge consists of reaching low enough tem- 
peratures/entropies for the observation of ordered mag- 
netic states. Theoretical insight on optimal conditions 
for such observations is greatly needed. While Mott in- 
sulator (MI) phases of single component bosonic systems 
have been observed experimentally 0, 0, [lfl | , and finite 
temperature effects have been extensively investigated re- 
cently [ill, [H, [TH, Q [H[ , the multi-component case is 
still a work in progress. 

In the present work, we address the issue for the case 
of two-component bosonic systems. We obtain such im- 
portant numbers as critical temperatures and, more im- 



portantly entropies, below which magnetic phases can be 
observed experimentally. With these numbers in hand, 
we provide rough estimates of hold times required for ob- 
serving thermally equilibrated ordered magnetic states. 

We consider a homogeneous system of two-component 
bosons in a cubic (square) lattice with repulsive inter- 
species interaction and half-integer filling of each compo- 
nent. This system can be realized by loading OL with 
two different atomic species, see, e.g., experiments at 
LENS with rubidium and potassium mixtures (T5 , [l7| , or 
the same atomic species in two different internal energy 
states, see e.g recent experiments done at MIT [l5| and 
ongoing experiments at Stony Brook [Lsj . The inter- and 
intra-species interaction strengths, U a b = U, U aa , and 
Ubb can be tuned via Feshbach resonance or by changing 
the Wannier functions overlap (in the presence of state- 
dependent lattices). If the intra-species interactions U aa 
and Ubb are made much larger than any other energy 
scale, and the temperature is low enough, the system 
is accurately described by the two-component hard-core 
Bose-Hubbard Hamiltonian: 
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Here aj^aj, b[{b^j are bosonic creation (annihilation) op- 
erators and t a , tb are hopping matrix elements for two 
species of bosons (A and B), respectively; the symbol 
< . . . > imposes the nearest-neibor constraint on the sum- 
mation over site subscripts; nf^ = a|a i and nf 1 ' = 6^6 i . 

Model (TIT) displays a very rich ground state phase di- 
agram [nj, [20, [21(, see Fig. [TJ For strong enough in- 
teractions, the system is incompressible in the particle- 
number sector, i.e. it is a MI. The remaining degree of 
freedom describing the boson type on a given site can be 
mapped onto the effective iso-spin variable [l^, H3, [HJ 
and gives rise to two possible MI states: a double checker- 
board (2CB) solid phase, equivalent to the Ising antifcr- 
romagnet, and a super-counter-fluid (SCF), equivalent 
to a planar ferromagnct in the iso-spin terminology. For 
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large enough hoppings the MI state undergoes a transi- 
tion to a double superfluid state (2SF). Finally, as it has 
been shown recently [2l[ , for strong asymmetry between 
the hopping amplitudes and relatively weak inter-species 
interaction a solid phase in the (heavy) component is 
stabilized via a mechanism of inter-site effective interac- 
tions mediated by the (light) superfluid component. In 
what follows we will focus on the magnetic states, namely 
the Ising antiferromagnet and the a;?/-fcrromagnct. We 
present the first precise results, based on path integral 
Monte Carlo (PIMC) simulations by the Worm algorithm 
[23l | . for transition lines to magnetic phases in two- and 
three-dimensions (2D and 3D) at zero and finite temper- 
ature, and discuss experimental parameters required for 
reaching them. 



II. GROUND STATE 

We begin with results for the ground state. In Fig. [T] 
we show the complete zero temperature phase diagram of 
model (H|) for the 2D system calculated in Ref. [2l|. We 
also sketch (dashed line) the transition line for the disap- 
pearance of magnetic order for the 3D system by comput- 
ing benchmark transition points (down triangles) for the 
strongly anisotropic and isotropic limits. These points 
correspond to the disappearance of the insulating Ising 
and the (ij/)-ferromagnetic phases, respectively. While, 
as expected, the 3D case is better captured by the mean- 
field theory [2(| HH , the discrepancy between mean-field 
and Monte Carlo results is still sizable: ^50%. 

These results provide quantitative guidance for experi- 
mentally achieving the regime of quantum magnetism. In 
experiments with two different species this can be done 
by using Feshbach resonances [l6| in order to reach the 
desired t a ,i,/U value; in the case of the same species but 
different internal states one can load state dependent lat- 
tices and tune the interspecies interaction by changing 
the overlap of Wannicr functions of the two components. 



III. FINITE-TEMPERATURE RESULTS 

Turning to the issue of reaching magnetic phases in re- 
alistic experimental setups — with an adiabatic protocol 
of turning on the optical lattice — we look for highest pos- 
sible values of the critical entropy for the appearance of 
magnetically ordered states. The critical values of tem- 
perature come as a natural 'by-product' of simulations. 
In what follows we use % > t a as the energy unit. 



A. Critical temperatures 

We start with the Ising antiferromagnet-to-normal 
transition. It belongs to the d-dimcnsional Ising univer- 
sality class, the order parameter being the staggered mag- 
netization along the z-axis or, equivalcntly, in bosonic 
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FIG. 1: (Colors online). Phase diagram of model ([TJ on a 
square lattice and half integer filling factor of each compo- 
nent (z is the coordination number). The 2CB-SCF first-order 
transition is represented by circles, the SCF-2SF second or- 
der transition by squares, the 2CB-2SF first-order transition 
by stars, the 2CB-(CB+SF) second-order transition by di- 
amonds, and the (CB+SF)-2SF first-order transition by up 
triangles. Down triangles are benchmark points for the dis- 
appearance of magnetic order in the cubic lattice. Lines are 
to guide an eye. 



language, the structure factor (which is the square of the 
order parameter): 



? (a,b) 
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with K the reciprocal lattice vector of the CB solid, i.e. 
K=(7r, 7r) in 2D and K=(7r, 7t,7t) in 3D, nf"^ the filling 
factor at the site r, and N^ a ' b ' 1 the total number of par- 
ticles A, B. In the vicinity of the transition point, the 
structure factor scales as 

2/3 2/3 1 

SK(T,L)=C-m/L) = L--g( T L»), (3) 

where £ is the correlation length, r = (T — T c )/tb is the 
reduced temperature, L is the system size, assumed to 
be large enough to neglect higher-order corrections to 
the universal scaling, f(x) and g{x) are universal scal- 
ing functions, and (3 and v are the critical exponent 
for the order parameter and correlation length, respec- 
tively. For the 2D case 2(3 jv = 1/4, and for the 3D case 
2(3 1 v = 1.0366(8) [Hj]. At the critical point, the quantity 
Sy^L 2 ^I v is size independent, provided L is appropriately 
large, and curves of different L's intersect. Figure[2]shows 
an example of the intersection for the case of a 2D sys- 
tem, with parameters t a /tb = 0.285 and U/% = 5.7, and 
system sizes L = 8, 16, 20, 24, 30. The critical tempera- 
ture is T c /t b = 0.1175(10). 

Our results for critical temperatures in 2D are sum- 
marized in Fig. [3] We have performed simulations at 
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FIG. 2: (Colors online). Finite size scaling for the struc- 
ture factor in the 2D system (see text) for t a /tb = 0.285, 
U/t b = 5.7 and system sizes L = 8,16,20,24,30. The criti- 
cal temperature can be read form the intersection of curves 
corresponding to different L's. Lines are a guide to an eye. 
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FIG. 3: (Colors online). Critical temperature for the Ising 
state vs. t a /U in the 2D system at fixed 2ztt/U = 0.1 and 
0.2, squares and circles respectively. Lines are a guide to an 
eye. 



fixed 2zt a /U = 0.1,0.2 and varying t^/U. Our data 
show that the region with higher transition temperatures 
corresponds to relatively weak interactions, but away 
from the transition to the (CB+SF) ground state. For 
strong interactions, the relevant energies, i.e. coupling of 
spin degrees of freedom in the mapping to the quantum 
spin Hamiltonian, scale as cx U -1 , and therefore require 
smaller temperatures in order to stabilize magnetically 
ordered phases. On the other hand, for weak enough in- 
teractions, the magnetic order will eventually disappear 
in favor of the (CB+SF) phase. As we approach this 
transition the magnetic order becomes weaker, therefore 
lower temperatures are required to observe it though the 
effect is rather moderate. The largest transition temper- 
atures lie somewhere in between this two limits, and with 
precise numerical simulations it is possible to accurately 
pinpoint the parameter region which is best suited for 
current experiments. The largest critical temperatures 
we have observed are T c /i& ~ 0.12. 

In the 3D case, we have calculated T c in the region 
where we expect it to be large, U/tb = 11, t a /tb = 0.1. 
We have found T c /t h = 0.175(15). The 3D simulations 
are far more demanding computationally than in 2D, and 
the calculation of the full zero- and finite-temperature 
phase diagram in 3D is beyond the scope of this work. 

We now turn to the melting of the xy-fcrromagnctic 
state. In bosonic language, it corresponds to the SCF- 
to-normal transition where SCF is characterized as the 
superfluid state with the composite order parameter de- 
scribing the condensate of pairs consisting of particles 
of one component and holes of the other one, with zero 
net particle flux. The transition is of the d-dimensional 
U(l) universality class, meaning that in 2D it is of the 
Kostcrlitz-Thouless (KT) type. In Fig. 2] we show an 
example of how transition points for the 2D system 




FIG. 4: (Colors online). Main plot: superfluid stiffness of the 
particle-hole composite object in the SCF (zy-ferromagnetic) 
state with t a = ii> and U/tb = 11 for system sizes L=10, 20,40 
trangles, circles, squares respectively. Inset: scaling of the 
finite-size 'critical temperature' (see text). 

are calculated. In order to locate the critical temper- 
ature we employ finite-size arguments following from KT 
renormalization-group flow for the superfluid stiffness p Sl 
the latter being measured from statistics of fluctuations 
of winding numbers [25j : 

» = W^> (4) 

where W is the vector of worldline winding numbers in 
the SCF sector. For our purposes, it is sufficient to define 
p s up to a global pre-factor; that is why our Eq. (H} 
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contains no other factors. 

In terms of worldlinc windings, the universal Nelson- 
Kostcrlitz jump translates into the abrupt change of 
(W 2 ) at the critical point from A/ir in the SCF phase 
to zero in the normal phase. In a finite system, the uni- 
versal jump is smoothed out and winding numbers go 
to zero continuously (see the main plot in Fig. 2]). If 
one defines the finite-size critical point T C (L) by the con- 
dition (W 2 (T C (L))) = 4/tt, then the flow of T C {L) to 
the thermodynamic limit answer T c = T c (oo) is given by 
T C (L) -T c (x l/(lni) 2 , see the inset in Fig. [2 

We have found the following critical temperatures: 
T c /t b = 0.141(5) for U/h = 11, t a /t b = 1; T c /t b = 
0.104(5) for U/t b = 13, t a /t b = 1; T c /t b = 0.101(5) for 
U/t b = 11, t a /h = 0.8; T c /t b = 0.14(1) for U/t b = 9.4, 
t a /t b = 0.6. Critical temperatures seem to decrease as 
we go towards the Heiscnberg point and the effective iso- 
spin couplings decrease (see argument above). Unlike 
the Ising-normal transition, the highest transition tem- 
perature we have found lies close to the SCF-2SF T = 
transition line. In fact, across this transition line the su- 
perfluid stiffness of the particle-hole composites and the 
transition temperature to the normal state remain finite. 
As discussed in Ref. [2f|, at finite temperature the SCF- 
2SF boundary moves in the direction of the 2SF ground 
state thus implying the following sequence of events: as 
temperature is increased in the vicinity of the quantum 
critical point the 2SF state first undergoes a transition to 
the SCF state which then turns normal at a much higher 
temperature. 

In the 3D case, the transition point can be obtained 
from the finite-size scaling of p s . Similarly to Eq. ([3]), 
one has: 

Ps (T,L)=C 1 f((/L)=L- 1 g(TLi). (5) 

The critical temperature is extracted from the intersec- 
tion of p s (r,L)L curves. We have done simulations for 
the system parameters U/t b = 21, t a /t b = 1 and found 
T c /t b = 0.208(7). 



B. Entropy curves 

Entropy curves S(T) are calculated starting from the 
energy data. We first use spline interpolation of data 
points to obtain a smooth curve E(T). We then calcu- 
late entropy by using two different numerical procedures: 
(i) We obtain the specific heat cy by differentiating the 
spline and then calculate the entropy by numerical inte- 
gration of cy/T. (ii) We avoid numerical derivatives by 
using 

and numerical integration. The agreement of the two 
methods is very good (within 0.5%). Uncertainties in 
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FIG. 5: (Colors online). Entropy curves for the Ising anti- 
ferromagnet in 2D for U/tb = 5.7, t a /tb = 0.142 and 3D for 
U/tb = 11, ta/tb — 0.1, solid and dashed lines respectively. 
Dotted lines are a guide to the reading of critical entropies. 




FIG. 6: (Colors online). Entropy curves for the xy- 
ferromagnet in 2D for U/tb = 11, t a /tb = 1 and 3D for 
U/tb = 21, ta/tb — 1, solid and dashed lines respectively. 
Dotted lines are a guide to the reading of critical entropies. 



entropies come therefore from the ones in critical temper- 
atures and finite-size effects. Examples of entropy curves 
in the Ising antifcrromagnctic state are shown in Fig. [5l 
for U/t b = 5.7, t a /t b = 0.1425 in 2D, and U/t b = 11, 
ta/h = 0.1 in 3D. We find critical entropies per parti- 
cle S c {k B )/N ~ 0.25 ± 5% and 0.5 ± 20% in 2D and 
3D, respectively. These entropies are relatively large and 
definitely within the realm of what can be achieved with 
bosonic BECs. In Fig. [5] we show entropy curves for 
the xy-ferromagnetic state. The critical entropy in 2D 
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for U/t b = 11, t a /t b = 1 is S c (k B )/N ~ 0.033 ± 5%, 
about an order of magnitude smaller (!!) than for the 3D 
value S c {k B )/N ~ 0.35 ± 10% obtained for U/t b = 21, 
ta/h = 1- This is explained by the specifics of the KT 
transition when the SF density jumps to zero discontin- 
uously at the critical point, i.e. when the system ther- 
modynamics is still dominated by the dilute phonon gas. 
Correspondingly, at the transition temperature the ther- 
mal energies and entropies are low. Our thermodynamic 
data confirm that this is precisely what is happening for 
the 2D system: energy scales with temperature as oc T 3 
(which implies that entropy is oc T 2 ) all the way up to 
temperatures T <T C . 



IV. MINIMAL EXPERIMENTAL HOLD TIMES 

Finally, we estimate minimal hold times required to 
observe ordered magnetic phases under typical experi- 
mental conditions. 

For a cubic lattice and using a harmonic approxima- 
tion around the minima of the optical lattice potential 
[27| . the tunnelling matrix elements and on-site interac- 
tion energies are given by: 
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co a , b = ^4E^ h) V a ' h /H 
is the harmonic oscillator frequency, 
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is the atomic recoil energy, m a , b and v ab are the bare 
and reduced masses respectively, a 'f a ' bb * > anc ] a < f b '> are 
the intra- and interspecies scattering lengths. The hard 
core limit can be achieved if e.g. ai aa " bb ^ < a[ ab \ or by 
manipulation of the overlapping of Wannicr functions as 
explained above. For 87 Rb - 41 K mixtures [HI and away 
from resonances one has aptb-K = 163ao, <XRb = 99ao, 
and ax = 65ao (do is the Bohr radius). One can then 
use Feshbach resonances to tune scattering lengths to the 
hard-core limit. 



To estimate the hold time t cxp required for the ob- 
servation of the magnetic phases we look at the lowest 
dynamic energy scale in the system which is t a in our 
case. Clearly, unless a condition i cxp 3> h/t a is satisfied, 
one may not even discuss thermally equilibrated normal 
states, not to mention low temperature ordered ones. If 
T c is smaller than t a , we consider t cxp S> h/T c as the 
minimal requirement. As we have seen, the optimal ex- 
perimental parameters for both Ising and xy phases result 
in min(£ Q ,T c ) ~ 0.1t b , and in what follows we will use 
this energy scale for the estimate of the hold time. 

Let us consider laser beams with A =1064 nm and dis- 
cuss the mixtures of Rb atoms in states | 1,-1) and 

I 2,-2) for which a * a6) = 98 - 09a o- For thc 

melting of the Ising state we require U jt b ~ 11 and 
ta/tb ~ 0.1 which translates into the optical lattice 



depths Va/Eg 1 ~ 19.5 and V b /E%> ~ 9, and the fi- 
nal result t cxp ^> 0.2s. For thc melting of thc xy- 
ferromagnet we require U/t b ~ 21, t a ~ t b , or, in 

terms of the lattice depths, V a /E^ ] = V b /E^ ] ~ 12, 
which implies that t cxp ^> 0.035s. For the case of 
87 Rb - 41 K mixtures, the best-case scenario corresponds 
to a.R,b-K = 163ao and aRb,ax flRb— k achievd via 
Feshbach resonances for intraspecies collisions. We con- 
sider the b species to be 87 Rb. A similar analysis of 
the Ising antifcrromagnctic case leads to V a /E$ ~ 19.5, 
V b /E^ ~ 6 and t cxp 3> 0.08s. For thc xy-ferromagnetic 

3.6 and 



,(o) 



(6) 



11.5 and V b /E% 



case we have V a /E R 
t, 

scattering length to, e.g., aRb-K ~ 35ao, in order to 

4 a) 



cxp 3> 0.015s. If, instead, one tunes the interspecies 
catl 

achieve the hard-core limit, this implies V a /E^ ~ 26.2, 
V b /E^ ~ 10.6, t cxp 0.25s for thc Ising antiferromag- 



net, and V a /E^ ] ~ 17.3, V b /E%> ~ 13.6, t cxp > 0.05s 
for the xy-ferromagnet. 

From these estimates we conclude that observing or- 
dered magnetic phases will be experimentally challenging 
since the required time scales might have to exceed sec- 
onds (with some advantage for dealing with the 87 Rb 
- 41 K mixture). Increasing the sample stability and sup- 
pressing various heating mechanisms (three-body losses, 
background vacuum, spontaneous scattering of lattice 
photons, and technical noises such as beam alignment, 
intensity fluctuations, mechanical vibrations) has to be 
achieved. To appreciate the problem, we mention the 
heating rate (entropy per particle) of ~ lfcs/s observed 
recently in a typical experiment in the optical lattice pd| . 



V. CONCLUSION 

We have addressed numerically (by worm algorithm 
Monte Carlo simulations) the problem of magnetic or- 
dering in the two-component Bose-Hubbard model in the 
intraspecies hard-core limit, for 2D and 3D cases, at fi- 
nite temperature. The emphasis of thc study is on re- 
vealing the optimal parameters for (and analyzing the 
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feasibility of) experimentally achieving the transitions to 
Ising antiferromagnetic (a.k.a. checkerboard solid) and 
zy-ferromagnetic (a.k.a. super-counter-fluid) phases. We 
have identified the optimal experimental conditions, cor- 
responding to maximal critical entropy per particle, tem- 
peratures and entropies. On the basis of our data, we 
have estimated minimal experimental hold times required 
to reach equilibrium magnetic states. These times have 
to be on a scale of seconds which renders the experimen- 
tal observations of magnetic phases challenging and calls 
for increased control over heating sources. 

Our results — optimal Hamiltonian parameters with 



corresponding values of critical entropies, temperatures, 
and minimal hold times — can be directly used for 
guiding and benchmarking the on-going experiment on 
creating optical lattice emulators. 
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